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Abstract 

The  objectives  of  localized  density  modeling  is  presented  and  compared  to  existing  global  models  currently  in  use. 
Localized  models  are  primarily  concerned  with  the  density  profiles  observed  by  a  specific  spacecraft  in  a  slow 
varying  orbit  regime.  The  CHAMP  density  and  derived  temperature  data  are  presented  and  discussed  as  relating  to 
localized  modeling  and  truth  values.  Temperature  perturbations  due  to  geomagnetic  storms  are  discussed  and  the 
basis  of  a  new  model  that  applies  local  information  is  discussed.  Output  based  models  are  developed  based  on  the 
autoregressive  filters,  that  use  previous  density  measurements  to  make  predictions  of  future  density  values.  Future 
work  of  combining  simplified  local  density  models  with  perturbation  corrections  is  also  discussed. 


1.  Introduction 

Since  the  development  of  Luigi  G.  Jacchia’s  first  density  model  in  1970  (J70)  [1],  atmospheric  density  modeling  has 
steadily  focused  on  large  monolithic  codes  that  provide  global  density  coverage  [2]  [3].  The  most  recent 
instantiation  of  the  global  density  model  is  the  Jacchia-Bowman  2008  (JB2008)  model  developed  by  Bruce  Bowman 
of  the  Air  Force  Space  Command  [4],  As  the  models  have  evolved  and  improved,  their  complexity  has  grown  as 
well.  Where  the  J70  model  required  2  indices  and  various  time  averages  to  determine  density,  the  JB08  model 
requires  5  indices  to  determine  density  [4],  Due  to  computational  complexity,  the  number  of  real-time  inputs 
required,  and  limited  forecasting  abilities,  these  models  are  not  well  suited  for  onboard  satellite  orbit  propagation. 

In  contrast  to  the  global  models,  this  paper  proposes  the  development  of  a  density  prediction  tool  that  is  only 
concerned  with  the  trajectory  of  a  specific  satellite.  Since  the  orbital  parameters  of  most  low  Earth  orbiting  satellites 
remain  relatively  constant  in  the  short  term,  there  is  also  minimal  variation  in  the  density  profile  observed  by  the 
satellite.  Limiting  the  density  model  to  a  smaller  orbit  regime  will  also  increase  the  ability  to  forecast  the  density 
along  that  orbital  track.  As  a  first  step,  this  paper  evaluates  the  feasibility  of  using  a  localized  density  prediction 
algorithm  to  generate  the  density  profile  that  will  be  seen  by  satellite,  allowing  for  high-accuracy  orbit  propagation 
with  minimal  or  no  input  from  the  ground. 

The  algorithm  evaluated  in  this  paper  is  a  simple  autoregressive  model  that,  given  previously  measured  density 
values,  provides  predictions  on  the  upcoming  density  profile.  This  first  approach  requires  zero  information  about 
the  satellite’s  current  orbit,  but  does  require  an  onboard  method  for  determining  the  current,  local  density.  Though 
this  aspect  of  the  onboard  system  is  not  analyzed  here,  it  is  envisioned  that  this  current,  local  density  (or  equivalently 
drag  acceleration)  would  be  calculated  through  onboard  processing  of  GPS  or  accelerometer  data. 

Using  the  trajectory  of  the  CHAMP  satellite  as  a  test  case,  various  samples  of  CHAMP  density  data  in  the  “past” 
(i.e.  before  a  chosen  epoch  time)  will  be  input  to  the  model,  and  the  model  will  in  turn  predict  “future”  density 
values  (i.e.  beyond  the  chosen  epoch  time).  The  effectiveness  of  the  local  density  model  will  be  assessed  by 
comparing  its  predicted  density  values  to  “true”  densities  in  the  CHAMP  database  at  those  times.  Results  describe 
the  prediction  accuracy  of  the  filter  and  length  of  time  over  which  accuracy  is  maintained. 
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2.  CHAMP  Density  and  Temperature  Data 


Launched  on  July  15,  2000,  the  CHAllenging  Minisatellite  Payload  (CHAMP)  has  provided  a  wealth  of  data  from 
its  high  accuracy  onboard  accelerometers.  Of  particular  interest  to  the  astrodynamics  community  is  the  atmospheric 
neutral  density  values  derived  from  the  drag  accelerations  [5]  [6],  The  derived  density  for  May  20,  2001  (day  140) 
are  shown  in  Fig.  la.  Fig.  lb  shows  the  derived  density  for  the  first  4  hours  of  the  same  day  along  with  the  densities 
predicted  by  the  Jacchia-Bowman  2006  (JB2006)  [2]  and  the  High  Accuracy  Satellite  Drag  Model  (HASDM)  [3] 
global  models.  Note  that  while  the  JB2006  and  HASDM  models  are  fairly  consistent  with  each  other,  they  both 
over-predict  the  density  for  the  time  shown.  For  this  analysis  the  CHAMP  derived  density  is  considered  to  be  truth. 


CHAMP  Density:  2002  Day  140 
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Fig.  1.  (a)  CHAMP  derived  neutral  density  (logarithmic  scale)  for  2001  Day  140.  (b)  CHAMP  density  data 
compared  to  the  JB2006  and  HASDM  predicted  densities. 


Another  set  of  data  derived  from  the  CHAMP  data  is  the  atmospheric  temperature  variations  due  to  geomagnetic 
perturbations.  For  the  J70  model  and  many  of  its  derivative  models  this  is  the  AT  associated  with  the  geomagnetic 
indices  ap  or  kp,  and  is  the  primary  means  of  addressing  geomagnetic  storms  and  their  related  density  fluctuations 
[7].  This  data  is  formed  by  evaluating  the  J70  model  for  all  of  the  satellite  position  and  known  indices,  minus  the  ap 
value.  The  algorithm  then  determines  the  AT  value  (typically  a  function  of  ap)  that  would  result  in  the  CHAMP 
observed  density.  The  AT  values,  as  a  function  of  the  ‘observed’  ap  indices,  for  years  2001  through  2005  are  shown 
in  Fig.  2.  Note  that  the  data  has  a  large  spread  for  each  value  of  ap.  Similarly  derived  data  was  used  to  determine 
the  new  DST  based  geomagnetic  perturbation  model  employed  in  the  JB2008  model  [8]. 


2001-2005  Storm  Data 


Fig.  2.  J70  geomagnetic  temperature  perturbations  derived  from  the  CHAMP  density  data.  The  red  line  is  the  AT 

predicted  by  Jacchi’s  static  model. 


3.  Geomagnetic  Storm  Perturbations 


In  an  effort  to  evaluate  the  applicability  of  localized  density  modeling,  an  effort  was  first  made  to  better  understand 
the  density  perturbations  due  solely  to  geomagnetic  storms  using  the  CHAMP  derived  AT  data.  This  work  was  also 
being  performed  as  an  ap  based  backup  to  the  new  DST  models  in  the  JB2008  global  density  model.  In  his  paper 
entitled  “New  Static  Models  of  the  Thermosphere  and  Exosphere  with  Empirical  Temperature  Profiles,”  Jacchia 
outlined  his  empirical  relation  between  the  Kp  and  ap  indices  and  AT  to  be  the  following  [7]: 

ATg  =  28  °Kp  +  0.03°eKv  (1) 

ATg  =  1.0°ap  +  100°[1  -  e-°  °8ap]  (2) 

Jacchia’s  model  is  global  in  formulation  and  does  not  account  for  atmospheric  variations  resulting  from  latitudinal, 
temporal,  or  seasonal  variations.  The  models  also  possess  an  inherent  average  time  lag  between  the  variations  in 
geomagnetic  index  and  those  in  temperature  of  6.7  hours. 

After  multiple  evaluations  and  iterations  modeling  the  data  an  more  localized  modeling  approach  was  developed  to 
represent  the  observed  data.  To  do  this  a  new  modeling  parameter,  0,  was  developed  to  indicate  both  the  latitude  of 
the  spacecraft  as  well  as  whether  its  current  location  was  on  the  day  or  night  sides  of  the  earth  and  in  the  summer  or 
winter  hemisphere.  Due  to  the  limitations  in  the  data  no  binning  was  performed  about  the  local  solar  time  (LST)  or 
day  of  year  (DOY).  The  parameter  0  is  depicted  in  Fig.  3.  Fig.  3a  shows  the  earth,  projected  as  a  circle,  divided  into 
four  quadrants  representing  the  unions  of  day/night  and  summer/winter.  Fig.  3b  depicts  the  same  projected  earth 
with  the  equivalent  0  angle  (which  is  simply  a  circle  with  angle  gradations  from  0  to  360).  A  satellite  in  the  summer 
hemisphere,  on  the  night  side  of  the  earth,  and  located  at  approximately  60  degrees  latitude  (whether  positive  or 
negative)  would  have  a  0  value  of  120  degrees.  A  satellite  in  the  winter  hemisphere,  on  the  day  side  of  the  earth, 
and  located  at  approximately  ±60  degrees  would  have  a  0  value  of  300  degrees. 

Since  the  objective  of  the  model  being  developed  was  an  ap  based  model,  the  data  was  first  sorted  according  to  its 
corresponding  ap  value,  and  then  binned  according  to  its  calculated  value  for  0.  The  results  are  shown  in  Fig.  4. 
Fig.  4a  depicts  the  AT  data  in  a  polar  plot  with  each  ap  set  normalized  to  its  maximum  value.  Fig.  4b  depicts  the  AT 
data  vs  0  for  each  discrete  value  of  ap  (including  ap  =  0). 
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Fig.  3.  (a)  Earth  divided  into  day/night/summer/winter  quadrants,  (b)  Model  parameter  0  depicting  appropriate 

quadrant  and  latitude. 
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Fig.  4.  (a)  Polar  plot  of  normalized  AT  data  for  each  discrete  value  of  ap.  (b)  AT  data  for  each  discrete  value  of  ap 

vs  the  parameter  0. 


It  is  not  the  intent  of  this  paper  to  interpret  the  physical  aspects  of  this  data  (such  as  the  perturbations  due  to  seasonal 
conditions)  but  to  analyze  the  ability  to  apply  corrections  to  a  localized  density  model.  An  inspection  of  the  data  in 
Fig.  4b  indicates  that  temperature  perturbations  (and  equivalently  density  perturbations)  due  to  magnetic  storms  do 
show,  to  some  extent,  predictable  behavior  about  an  orbit.  The  data  for  each  value  of  ap  can  easily  be  modeled  by 
a  second-order  fourier  series.  Given  a  base  model  for  the  local  density  it  should  be  possible  to  apply  correctors, 
based  on  observed  data,  to  provide  reasonably  accurate  localized  density  predictions.  The  remainder  of  this  paper 
address  the  some  of  the  initial,  though  incomplete,  steps  towards  developing  the  localized  base  model. 

4.  Localized  Model:  Output  Based  Predictions 


As  a  first  step  towards  a  base  model  we  evaluated  output  based  predictions,  based  on  signal  processing  theory. 
Signal  processing  filters,  such  as  the  autoregressive  moving  average  (ARMA)  model,  can  make  signal  predictions 
based  on  past  inputs  and/or  past  outputs  [9],  The  autoregressive  (AR)  model  makes  predictions  of  future  signal 
values  based  solely  on  past  outputs.  Using  this  type  of  model  allows  for  density  predictions  based  solely  on  past 
density  measurements  (derived  through  onboard  accelerometers).  The  generic  AR  model  is: 


yt 


n 


i= 1 


+  et 


(3) 


where  c  is  a  constant,  et  is  white  noise,  at  are  the  model  coefficients,  yt_j  are  previously  measured  signals,  and  yt  is 
the  predicted  signal.  Solving  for  the  model  coefficients  can  be  performed  using  Matlab’s  LPC  function. 


This  AR  model  predicts  one  step  into  the  future.  Given  that  the  CHAMP  data  samples  being  analyzed  are  spaced  80 
seconds  apart,  this  represents  a  model  that  has  an  80  second  forward  looking  prediction  capability.  The  results, 
shown  in  Fig.  5,  depict  that  these  short  term  predictions  work  reasonably  well.  With  this  type  of  based  model  there 
is  no  need  to  account  for  density  perturbations  beyond  the  basic  AR  model.  Predictions  a  few  minutes  into  the 
future  though  is  perhaps  a  trivial  solution  and  not  particularly  useful  for  long  term  predictions. 


Fit  of  Density  to  2nd-order  AR  (2002  Day  140) 


Fig.  5.  Fit  of  the  CHAMP  data  using  a  second-order  AR  model  to  predict  one  step  into  the  future.  Data  samples  are 

at  80  second  intervals. 


5.  Multi-Step  Predictions 

Now  that  we’ve  shown  the  ability  to  predict  the  density  one  step  into  the  future,  let’s  suppose  that  I’m  interested  in 
predicting  the  density  50  steps  into  the  future.  Reformulating  equation  3,  for  a  second-order  filter  yields: 

yt+so  =  brit  +  +  a2yt-2  (4) 

yt  =  brjt  +  +  a2yt_52  (5) 

This  indicates  that  data  being  predicted  is  the  linear  sum  of  data  51  and  52  time  steps  in  the  past.  Since  Matlab’s 
LPC  and  FILTER  functions  operates  on  the  idea  of  equally  spaced  data,  another  method  must  be  used  to  generate 
the  model  coefficients.  To  generate  the  coefficients  we  rewrite  the  equation  as: 

yt  -  [Tt-51  yt-52]  [a*]  +  brjt  (6) 

This  can  be  formulated  as  a  least-squares  problem: 

y  —  Hx  +  e  (7) 

Where  y  —  yt,  H  —  ^j,  x  —  | ai  a2 1,  and  e  —  br]t.  Solving  for  the  coefficients  x  provides  the  following 

equation: 


x  =  (. HTH)-1HTy  (8) 

Using  the  CHAMP  density  data  from  May  20,  2001  (DOY  140)  to  predict  73  time  steps  into  the  future 
(approximately  one  orbital  period  at  80  second  data  samples)  yields  the  following  values  for  the  filter  coefficients. 


ax  =  0.19453 
a2  —  0.80313 


Given  one  orbit  worth  of  density  data  (approximately  72  values  at  the  given  data  rate)  we  can  now  predict  the 
density  over  the  next  orbit.  The  model  is  then  applied  to  the  following  days  worth  of  data,  using  the  first  72  values 
to  seed  the  model.  The  results  for  May  21,  2001  (DOY  141)  are  shown  in  Fig.  6.  Note  that  the  future  predicted 
orbit  densities  are  similar  (though  not  identical)  to  the  previously  measured  density  values.  Inspection  of  the  results 
shown  that  this  model  responds  to  density  perturbations  with  a  one  orbit  delay.  The  model  is  essentially  saying  that 
without  any  additional  information  the  density  profile  of  the  next  orbit  should  look  similar  to  the  previous  density 
profile.  The  advantage  of  this  model  is  that  we  can  model  (at  a  base  level)  the  density  profile  for  one  full  orbit  as 
opposed  to  a  few  minutes  for  the  previous  formulation.  This  same  formulation  can  also  be  used  to  predict  multiple 
orbits  into  the  future  by  simply  repeating  the  predicted  orbit  profile  multiple  times. 

Fig.  7  depicts  the  ratios  of  the  CHAMP  data  to  the  JB2006,  HASDM,  and  one  orbit  AR  models.  Note  that  the 
JB2006  and  HASDM  models  both  over-predicted  the  density,  while  the  AR  model  has  a  mean  ratio  near  zero.  The 
JB2006  and  HASDM  models  are  both  operating  in  an  open  loop  while  the  AR  model  uses  a  moving  window  of  the 
‘measured’  density  to  model  future  values. 


Single  Orbit  Density  Prediction  (2001  Day  141) 


Fig.  6.  Fit  of  the  CHAMP  data  using  a  second-order  AR  model  to  predict  one  orbit  into  the  future.  Data  samples 

are  at  80  second  intervals. 
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Fig.  7.  Ratios  of  the  CHAMP  density  data  to  the  JB2006,  HASDM,  and  AR  model  predictions. 


It  can  be  pointed  out  that  the  results  shown  in  Fig.  6  are  for  fairly  benign  geomagnetic  conditions  and  are  based  on 
coefficients  generated  on  the  previous  days  data.  To  determine  if  this  model  approach  can  apply  to  perturbed 
condition,  such  as  during  a  geomagnetic  storm,  the  model  was  re  applied  to  October  3,  2001  (DOY  276).  This  day 
experienced  a  strong  geomagnetic  storm  with  an  ap  value  of  132  (Kp  =  7).  The  same  coefficients  are  maintained  for 
this  new  data  set.  The  AR  model  results  are  shown  in  Fig.  8.  The  ratios  for  the  JB2006,  HASDM,  and  AR  models 
are  shown  in  Fig.  9. 
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Fig.  8.  Fit  of  the  CHAMP  data  using  a  second-order  AR  model  to  predict  one  orbit  into  the  future.  Data  samples 

are  at  80  second  intervals. 
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Fig.  9.  Ratios  of  the  CHAMP  density  data  to  the  JB2006,  HASDM,  and  AR  model  predictions. 


The  AR  model  results  for  the  stormy  DOY  276  perform  similarly  to  those  seen  for  the  quiet  DOY  141.  The 
predicted  density  profile  for  the  following  orbit  are  similar  to  those  measured  in  the  previous  orbit.  Degraded 
performance  is  seen  as  a  result  of  the  significant  perturbations.  The  mean  of  the  predictions  for  DOY  141  and  276 
are  both  nearly  1,  but  DOY276  has  a  standard  deviation  of  13.8%  compared  to  8.2%  for  DOY  141. 


6.  Conclusion 

In  contrast  to  large  global  density  models,  a  new  approach  to  onboard  localized  density  modeling  is  sought. 
Localized  density  models  and  predictors  can  be  used  for  improved  onboard  propagators  where  drag  plays  a 
significant  role  in  A  new  approach  to  temperature  perturbations  resulting  from  geomagnetic  storms  was  presented. 
Results,  based  on  temperature  variations  derived  from  CHAMP  data,  show  that  storm  perturbations  do  show 
structure  and  predictability  and  can  likely  be  applied  to  base  level  models  to  account  for  variations  in  a  spacecraft’s 
local  density  profile.  An  approach  for  a  base  level  density  model  was  then  shown  which  uses  previously  measured 
density  values  in  an  autoregressive  model  to  predict  future  density  models.  The  autoregressive  model  was  able  to 
accurately  predict  density  values  in  the  short  term  (or  very  small  fractions  of  an  orbit).  The  model  was  modified  to 
predict  density  values  over  one  orbit  into  the  future  based  on  the  previous  orbit’s  density.  Results  showed  that  the 
model  predicts  the  future  orbit  density  profile  to  look  similar  to  the  previously  measured  profile.  The  model  also 
lagged  by  one  orbit  in  responding  to  density  perturbations  for  small  to  moderate  perturbations.  The  model  was  then 
run  against  data  during  a  strong  geomagnetic  storm.  While  the  model  behaved  similarly  to  the  case  with  small 
perturbations,  predicting  the  next  orbit’s  density  profile  will  be  similarly  to  the  previously  measured  profile,  the 
strong,  short-term  perturbations  in  density  were  not  accurately  captured.  Future  efforts  will  included  improved  base 
models  for  a  spacecraft’s  local  density  and  combining  perturbation  corrections  to  the  base  models.  Another  effort 
for  study  is  switching  from  a  local  density  model  to  a  local  drag  model,  eliminating  the  need  to  determine  density 
values  from  accelerometers  and  allowing  for  orientation  specific  accommodations  to  drag. 
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